function sandbox = fsb_get_hemodynamics_blocked(d,sandbox)

% FSB: Block design selection GUI
%
% EXAMPLE: sb = fsb_get_hemodynamics_blocked(150)
%
% INPUT:
% d: number of volumes
% TR: TR in seconds
%
% OUTPUT:
% sb: sandbox struct for blocked design
%
% CALLED BY:
% fsb_get_hemodynamics.m
%
% NOTES:
% For very basic on/off designs 
%
% Copyright 2010 MPI for Biological Cybernetics
% Author: Steffen Stoewer
% License:GNU GPL, no express or implied warranties
% 
% $ Revision 1.0
%
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

%~~~~~~~~~~~~~~~~~~~~~~~~~~
% Check input arguments
%~~~~~~~~~~~~~~~~~~~~~~~~~~
if nargin<2
    TR = 2;
    disp('No TR data provided, assuming TR = 2s');
else
    TR = sandbox.TR;
    disp(['TR is ' num2str(TR) ' s']);
end

%~~~~~~~~~~~~~~~~~~~~~~~~~~
% Set GUI up
%~~~~~~~~~~~~~~~~~~~~~~~~~~

selwindow = figure('Visible','off','Position',[100,100,480,260],'Resize','off','Toolbar','none','MenuBar','none','NumberTitle','off','Name','Blocked design setup');
figure(selwindow);

%~~~~~~~~~~~~~~~~~~~~~~~~~~
% Define blockdata struct
%~~~~~~~~~~~~~~~~~~~~~~~~~~

blockdata.blockon = 8;
blockdata.blockoff = 8;
blockdata.blockselect = 0;
blockdata.blockrep = 4;
blockdata.dispderiv = 0;

%~~~~~~~~~~~~~~~~~~~~~~~~~~
% Add Buttons and edit fields
%~~~~~~~~~~~~~~~~~~~~~~~~~~
uicontrol('Style','Text',...
    'Position',[300,230,50,25],...
    'String','TR:');


h_TR = uicontrol('Style','Edit',...
    'Position',[340,230,50,25],...
    'String',num2str(TR),...
    'Callback',{@TR_Callback});

uicontrol('Style','Pushbutton',...
    'Position',[10,140,100,25],...
    'String','Done',...
    'TooltipString','Done with setup',...
    'Callback',{@blockdone_Callback});

uicontrol('Style','Pushbutton',...
    'Position',[10,170,100,25],...
    'String','Preview',...
    'TooltipString','Show model',...
    'Callback',{@showmodel_Callback});

uicontrol ('Style','Text',...
    'Position',[150,170,75,25],...
    'String','ON per block',...
    'TooltipString','Number of volumes on per block')

h_blockon = uicontrol('Style','edit',...
    'Position',[230,170,25,25],...
    'String','8',...
    'TooltipString','Number of volumes on per block',...
    'Callback',{@blockon_Callback});

uicontrol ('Style','Text',...
    'Position',[150,200,75,25],...
    'String','OFF per block',...
    'TooltipString','Number of volumes off per block')

h_blockoff = uicontrol('Style','edit',...
    'Position',[230,200,25,25],...
    'String','8',...
    'TooltipString','Number of volumes off per block',...
    'Callback',{@blockoff_Callback});

uicontrol ('Style','Text',...
    'Position',[150,230,75,25],...
    'String','Repetitions',...
    'TooltipString','Number of block repetitions')

h_blockrep = uicontrol('Style','edit',...
    'Position',[230,230,25,25],...
    'String','4',...
    'TooltipString','Number of block repetitions',...
    'Callback',{@blockrep_Callback});

h_blockselect = uicontrol('Style','Radiobutton',...
    'Position',[10,200,100,25],...
    'String','Start with on',...
    'TooltipString','Start with on or off',...
    'Callback',{@blockselect_Callback});

h_dispderiv = uicontrol('Style','Radiobutton',...
    'Position',[10,230,100,25],...
    'String','Derivatives',...
    'TooltipString','Time and dispersion derivatives',...
    'Callback',{@dispderiv_Callback});

%~~~~~~~~~~~~~~~~~~~~~~~~~~
% Add display axis
%~~~~~~~~~~~~~~~~~~~~~~~~~~

graphax = axes('Units','Pixels','Position',[20,20,450,100]);

%~~~~~~~~~~~~~~~~~~~~~~~~~~
% Wait until user input is done
%~~~~~~~~~~~~~~~~~~~~~~~~~~

uiwait;

% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Button callbacks
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% TR callback
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

    function TR_Callback(source, eventdata)
        TR = get(h_TR,'String');
        set(h_TR,'String',num2str(TR));
        sandbox.TR = str2num(TR);
        TR = sandbox.TR;
    end

% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Block onset callback
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    function blockon_Callback(source, eventdata)
        blockon = get(h_blockon,'String');
        set(h_blockon,'String',num2str(blockon));
        blockdata.blockon = str2double(blockon);
                check_plausibility;
    end

% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Block offset callback
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    function blockoff_Callback(source, eventdata)
        blockoff = get(h_blockoff,'String');
        set(h_blockoff,'String',num2str(blockoff));
        blockdata.blockoff = str2double(blockoff);
                check_plausibility;
    end

% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Start with on or off callback
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

    function blockselect_Callback(source, eventdata)

        if (get(h_blockselect,'Value') == get(h_blockselect,'Max'))
            blockdata.blockselect = 1;
        else
            blockdata.blockselect = 0;
        end

    end

% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Show time and dispersion derivatives callback
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    function dispderiv_Callback(source, eventdata)

        if (get(h_dispderiv,'Value') == get(h_dispderiv,'Max'))
            blockdata.dispderiv = 1;
        else
            blockdata.dispderiv = 0;
        end

    end

% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Number of block repetitions callback
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    function blockrep_Callback(source, eventdata)
        blockrep = get(h_blockrep,'String');
        set(h_blockrep,'String',num2str(blockrep));
        blockdata.blockrep = str2double(blockrep);
                check_plausibility;
    end

% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Release lock and return callback
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    function blockdone_Callback(source, eventdata)
        uiresume;
        close(selwindow);
    end

% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Check plausibility
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    function check_plausibility(source, eventdata)
        if blockdata.blockrep*(blockdata.blockon+blockdata.blockoff)>d
            disp('length of experiment exceeds length of scan, aborting...')
                        set(gcf,'CurrentAxes',graphax);cla;
            text (0.2,0.5,'Design size does not match scan size, aborting');
            return
        end
    end

% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Preview model callback
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

    function showmodel_Callback(source,eventdata)
        
        try
        if blockdata.blockselect == 1;
            evt_onsets{1} = (blockdata.blockon + blockdata.blockoff).* [1:blockdata.blockrep] - (blockdata.blockon + blockdata.blockoff) + 1;
        else
            evt_onsets{1} = (blockdata.blockon + blockdata.blockoff).* [1:blockdata.blockrep] - blockdata.blockon +1;
        end

        %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
        % define Event onsets and event names
        %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

        evt_names{1} = 'stimon';
        evt_dur(1) = blockdata.blockon;

        onset = zeros(1,d);
        onset(evt_onsets{1}) = 1;

        evt_durreal{1} = evt_dur(1);

        inparg.evt_onsets = evt_onsets;
        inparg.evt_names = evt_names;
        inparg.evt_dur = evt_dur;
        inparg.TR = TR;
        inparg.volumes = d;
        inparg.evt_durreal = evt_durreal;
        inparg.deriv = blockdata.dispderiv;

        [U,bf,p,X] = fsb_conv_onsets(inparg);
        sandbox.expinfo = [];
        sandbox.hemodynamics = X;
        sandbox.onset = onset';
        sandbox.volmat.evt_onsets = evt_onsets;
        sandbox.volmat.evt_dur = evt_dur;
        
        %~~~~~~~~~~~~~~~~~~~~~~~~~
        % Setup intrial array
        %~~~~~~~~~~~~~~~~~~~~~~~~~
        sandbox.intrial = [];
        
        
            sandbox.intrial(:,1) = ones(d,1);
            sandbox.intrial(:,2) = repmat([1:blockdata.blockon+blockdata.blockoff],1,blockdata.blockrep);
            sandbox.intrial(:,3) = 1;
        catch
            set(gcf,'CurrentAxes',graphax);cla;
            text (0.1,0.5,'Design size does not match scan size, aborting');
            disp('Design size does not match scan size, aborting')
            return
        end

        y = 1;
        for x = 1:length(sandbox.intrial)-1
            if sandbox.intrial(x,2)>sandbox.intrial(x+1,2);
                y = y+1;
            end
            sandbox.intrial(x,3) = y;
        end
        sandbox.intrial(x+1,3) = y;

        %~~~~~~~~~~~~~~~~~~~~~~~~~
        % Keep everything compatible with FSB
        %~~~~~~~~~~~~~~~~~~~~~~~~~
        sandbox.timemat.evt_onsets = {};
        evt_onset = cellfun(@conv_time,evt_onsets,'UniformOutput',0);
        sandbox.timemat.evt_onsets = evt_onset;
       % evt_durs = cellfun (@conv_time,evt_dur,'UniformOutput',0);
        sandbox.timemat.evt_dur = ones(1,size(evt_onset{1},2))*evt_dur;
        sandbox.otto = [];
        sandbox.discard = [];
        sandbox.origtrial = max(sandbox.intrial(:,3));
        %sandbox.sensdata = [];
        set(gcf,'CurrentAxes',graphax);cla;
        bar(onset);
        hold on;
        plot(X,'r');
        axis tight;
    end

    function y = conv_time(x)
        y = x.*TR*1000;
    end

end
